
clear all

disp('Start');

addpath ../lib

fontSize = 26;

spektras1=load('c5.dat');
spektras2=load('fluorescein.dat');
spektras3=load('merocyanin.dat');
lamda=450:800;

%interpoliuojam, kad butu vienodi bangos ilgiai ir sunormuojam
sp1=interp1(spektras1(:,1), spektras1(:,2), lamda);
sp1(isnan(sp1)) = 0;
sp1=sp1/max(sp1);

sp2=interp1(spektras2(:,1), spektras2(:,2), lamda);
sp2(isnan(sp2)) = 0;
sp2=sp2/max(sp2);

sp3=interp1(spektras3(:,1), spektras3(:,2), lamda);
sp3(isnan(sp3)) = 0;
sp3=sp3/max(sp3);

h = figure('visible','off');
set(h, 'DefaultLineLineWidth', 1.5);
set(h, 'DefaultAxesFontSize', fontSize);

plot(lamda, sp1, '-b', lamda, sp2, '-r', lamda, sp3, '-g'); grid on;

xlabel('Bangos ilgis, nm', 'FontSize', fontSize);
ylabel('Intensyvumas, s.v.', 'FontSize', fontSize);
SaveFigure(h, 'markerSpectre');


sigma=5;
nwav=length(lamda);
nfilt=10;
npix=100;
filtras=zeros(10,length(lamda));
for i1=1:nfilt;
    lamdac(i1)=550+12*i1;
    filtras(i1,:)=exp(-((lamda-lamdac(i1)).^2)/(2*sigma^2));
end;

h = figure('visible','off');
set(h, 'DefaultLineLineWidth', 1.5);
set(h, 'DefaultAxesFontSize', fontSize);
plot(lamda, filtras); grid on;
xlabel('Bangos ilgis, nm', 'FontSize', fontSize);
ylabel('K, s.v.', 'FontSize', fontSize);
SaveFigure(h, 'filterSpectre');


image1=double(imread('rose.bmp')); image1=squeeze(image1(:,:,1));
image2=double(imread('tulip.bmp')); image2=squeeze(image2(:,:,1));
image3=double(imread('herbera.bmp')); image3=squeeze(image3(:,:,1));


%% -----
h = figure('visible','off');
subplot(2,3, 1);
imshow(image1,[], 'Border', 'tight');

subplot(2,3, 2);
imshow(image2,[], 'Border', 'tight');

subplot(2,3, 3);
imshow(image3,[], 'Border', 'tight');

%% ---

thres = 1;
image1 = Invert(image1);
image1 = image1 - min(min(image1));
image1(image1 < thres) = 0;

image2 = Invert(image2);
image2 = image2 - min(min(image2));
image2(image2 < thres) = 0;

image3 = Invert(image3);
image3 = image3 - min(min(image3));
image3(image3 < thres) = 0;


subplot(2,3, 4);
imshow(image1,[], 'Border', 'tight');

subplot(2,3, 5);
imshow(image2,[], 'Border', 'tight');

subplot(2,3, 6);
imshow(image3,[], 'Border', 'tight');

SaveFigure(h, 'initialPics');
%% ---

Data=zeros(nfilt,npix,npix);
%paveikslo inasas i-ojo filtro matavime - spektras x filtro pralaidumas ir 
%suintegruota
for i1=1:nfilt
    weight(i1,1)=sum(filtras(i1,:).*sp1);
    weight(i1,2)=sum(filtras(i1,:).*sp2);
    weight(i1,3)=sum(filtras(i1,:).*sp3);
    D = image1*weight(i1,1) + image2*weight(i1,2) + image3*weight(i1,3);
    Data(i1,:,:) = D;
end;




h = figure('visible','off');
for i1=1:nfilt
    subaxis(2, nfilt/2, i1, 'Spacing', 0.03, 'Padding', 0, 'Margin', 0);
    A=squeeze(Data(i1,:,:));
    imshow(Invert(A),[], 'Border', 'tight');
    title([int2str(lamdac(i1)) ' nm' ], 'FontSize', 20);
end;
SaveFigure(h, 'testiniaiDuomenys');


h = figure('visible','off');
for i1=1:nfilt
    subaxis(2, nfilt/2, i1, 'Spacing', 0.03, 'Padding', 0, 'Margin', 0);
    A=squeeze(Data(i1,:,:));
    imshow(A,[], 'Border', 'tight');
    title([int2str(lamdac(i1)) ' nm' ], 'FontSize', 20);
end;
SaveFigure(h, 'testiniaiDuomenysInvert');

h = figure('visible','off');
for i1=1:nfilt
    subaxis(2, nfilt/2, i1, 'Spacing', 0.03, 'Padding', 0, 'Margin', 0);
    A=squeeze(Data(i1,:,:));
    imshow(Invert(A),[], 'Border', 'tight');
end;
SaveFigure(h, 'testiniaiDuomenysNoTitle');


h = figure('visible','off');
for i1=1:nfilt
    subaxis(2, nfilt/2, i1, 'Spacing', 0.03, 'Padding', 0, 'Margin', 0);
    A=squeeze(Data(i1,:,:));
    imshow(A,[], 'Border', 'tight');
end;
SaveFigure(h, 'testiniaiDuomenysInvertNoTitle');

% Sukuriame galutine duomenu matrica RetData
% Kiekviena eilute tai isskleistas paveiksliukas vektoriumi su vis kitu
% filtru
for i1=1:nfilt
    A = squeeze(Data(i1,:,:));
    A = reshape(A,1,[]);
    RetData(i1,:) = A;
end


RetData = RetData';

%nustatome su SVD rango skaiciu
S = svd(RetData);
r = numel(S(S > 0.01));
%paleidziame Matlab'o nmf
m = max(max(RetData));
%RetData = RetData/m * 255;
%RetData = Invert(RetData);
[W, H] = nnmf(RetData, r);



%paplotiname inasus
x = 1:10;
h = figure('visible','off');
set(h, 'DefaultLineLineWidth', 1.5);
set(h, 'DefaultAxesFontSize', fontSize);
plot(x, H);
xlabel('Paveiksliuko Nr.', 'FontSize', fontSize);
ylabel('K', 'FontSize', fontSize);
SaveFigure(h, 'nmfWeights');

%plotiname isskirtus paveiksliukus
for i1=1:r
    C = W(:, i1);
    img = reshape(C, 100, 100);
    %img = 255 - img;
    img = Invert(img);
    
    h = figure('visible','off');
    imshow(img,[], 'Border', 'tight');
    SaveFigure(h, ['nmfPictures' int2str(i1)]);
end







disp('End');


    


    
